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ABSTRACT 

The  equations  obtained  by  Chien  for  the  non-linear 
deflection  of  shallow  spherical  shells  under  uniform 
external  pressure  are  solved  by  means  of  power  series  ex- 
pansions following  procedures  introduced  by  Frledrichs 
and  Stoker  in  th-lr  treatment  of  buckling  of  circular 
plates.   These  equations  depend  upon  tvD  parameters.   One 
of  these  is  related  to  the  external  pressure,  while  the 
other  depends  upon  the  dimensions  of  the  shell.   The 
equations  are  solved  for  several  ranges  of  the  parameters 
under  boundary  conditions  corresponding  to  a  fixed  edge. 

The  solution,  carried  out  numerically  on  the  Univac 
at  New  York  University,  yields  a  complete  description  of 
the  stresses  and  deflections  as  functions  of  the  polar 
angle  over  a  wide  range  of  values  of  the  loading  para- 
meter and  the  dimensional  parameter.   Prediction  of  the 
upper  buckling  load  is  then  made  by  means  of  a  numerical 
criterion  based  on  the  load  vs  deflection  curve.   For  some 
cases,  the  post-buckling  behavior  is  investigated.   The 
results  agree  well  with  existing  experimental  and  theo- 
retical studies  and  cover  a  wide  range  of  cases  not 
previously  treated. 
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SYMBOLS 

?EVi 
D       =  i-lexural   rigidity   of   shell     7- 

E  =  Young's  Modulus 

F,G  =  Boundary  functions,  see  Equation  19 

2h  =  Shell  thickness 

k  =  2/3(l-v^) 

N  =  Number  of  terms  In  power  series  required  for  convergence 

p  =  External  pressure,  positive  Inwards 

P  =  Loading  parameter,  see  Equation  12 

r  =  Radial  variable:   distance  to  any  point  in  shell 

R  =  Radius  of  middle  surface  of  shell 

S  =  Surface  of  shell 

u  =  Tangential  component  of  displacement 

u°  =  Tangential  component  of  displacement  on  middle  surface 

Y  =  Volume  of  shell 

w  =  Radial  component  of  displacement 

w°  =  Radial  component  of  displacement  on  middle  surface 

a  =  Slope  of  tangent  plane  to  deformed  surface  with  the 

horizontal  plane 

a  =  a/A 

Y  =  Stress  function,  see  Equation  10 

Y  =  Reduced  stress  function,  see  Equation  12 
e-  =  Longitudinal  strain 

e 

ex  =  ClrciAmTerential  strain 
e   =  A  small  number 
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e.  .  =  Components  of  strain  tensor 

K,  =  R-r 

\  =  A  small  number 

9  =  Polar  angle 

e  =  e/A 

=  Semi-angle  of  shell  opening 

X  =  Kaplan  and  Pung  geometric  parameter 

V  =  Poisson's  ratio 

17  =  Potential  energy 

p  =  Geometric  parameter,  see  Equation  12 

CTgjO^  =  Longitudinal  stress;  value  on  middle  surface 

^s»<yz  =  Circumferential  stress:  value  on  middle  surface 
9     9 

X  "^  =  components  of  stress  tensor 

^  =  Azimuthal  angle 
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Section  1»   Introduction 

Consider  a  thin,  hollow  sphere  of  wall  thickness 
2h  and  mean  radius   R  loaded  by  a  uniform  external 
pressure   p  directed  towards  the  center.   It  is  observed 
experimentally  that  as  the  press\ire  increases  from  zero, 
the  shell  apparently  deforms  by  contracting  into  a  smaller 
sphere.   This  mode  of  deformation  continues  vmtil  some 
critical  value  of  the  pressure  is  reached.   At  this  value 
the  shell  suddenly  jiamps,  or  buckles  into  a  non-spher- 
ical buckled  state  consisting  of  a  round,  rotatiorally 
symmetric,  inward  "dimple"  of  small  solid  angle  at  some 
point  on  the  svirface  of  the  shell,  while  the  remainder 
of  the  shell  remains  apparently  spherical. 

These  experimentally  observed  facts  are  entirely  at 
variance  with  the  linear  theory  of  buckling   .   Charac- 
teristically, the  linear  theory  predicts  a  wave  form  ex- 
tending over  the  whole  surface,  with  the  amplitude  unde- 
termined (e.g.  the  shell  could  buckle  inward  or  outward 
under  the  same  load).   This  is  to  be  contrasted  with  the 
experimentally  observed  "dimpling"  behavior.   Furthermore, 
the  experimental  buckling  load  is  found,  in  some  cases,  to 
be  as  low  as  one-quarter  of  that  predicted  by  the  linear 
buckling  theory. 

In  a  similar  type  of  phenomenon  observed  in  the  case 
of  buckling  of  thin  cylindrical  shells  under  axial  com- 
pression, attempts  were  made  to  explain  the  discrepancy 
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between  theory  and  experiment  on  the  basis  of  initial 
Imperf actions  from  the  cylindrical  shape  of  the  model, 
and  deviation  of  the  boundary  conditions  from  those  as- 
sumed in  the  theoretical  calculations.   Although  impor- 
tant, these  effects  were  not  sufficient  to  account  for 

the  lack  of  agreement  (see  the  discussion  In  reference  2). 

[2] 
Von  Karraan  and  Tslen    recognized  that  this  must  be 

the  result  of  a  non-linear  behavior  of  the  structure. 
Based  upon  the  minimization  of  a  potential  energy  ex- 
pression, they  shoivfed  that  for  values  of   p  belov/  the 
buckling  load  predicted  by  the  linear  theory  there  exist 
other,  non-adjacent,  stable  states  of  equilibrium.   The 
greatest  lower  bound  of  such  loads   p   is  called  the 
"minimum"  or  "lower"  buckling  load  p,..   If  p  exceeds 

Pt^  there  is  a  possibility  that  buckling  will  occur. 

[31 
Frledrlchs    pointed  out  some  shortcomings,  and 

suggested  improvements  in  the  Karman-Tsien  theory.   He 

also  treated  the  buckling  of  a  thin  sphere  as  a  "boundary 

layer"  phenomenon  where  the  "boundary  layer"  is  considered 

as  occuring  in  the  dimpled  region  of  the  sphere.   Tsien^  . 

adopting  some  of  the  suggestions  of  Frledrlchs  -^  ,  pre- 

r  pi 
sented  an  improved  version  of  the  Karman-Tslen    theory. 

[51 
Altshuler    attempted  to  determine  the  buckling  load  by 

solvinr  the  differential  equations  by  asymptotic  expan- 
sion techniques.  Although  none  of  these  approaches  are 
entirely  satisfactory,  they  do  reveal  the  essential  nature 
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of  the  phenomena  that  are  Involved  in  the  buckling 
process. 

Closely  r'elated  to  the  problem  of  the  buckling  of  a 
thin  sphere   is  that  of  the  buckling  of  a  shallow   spher- 
ical shell  under  uniform  external  pressure.   Physically 
as  well  as  mathematically   both  problems  reveal  a  simi- 
lar non-linear  behavior.   The  shallow  spherical  shoU  is 
of  interest  both  as  a  structural  entity  in  itself  and 
as  representing  that  portion  of  a  complete  sphere  in 
which  the  dimpling  occurs. 

Solutions  for  the  non-linear  deformations  of  shallow 

spherical  shells  under  various  external  pressures  and 

fip]  n^i 

edge  moments  have  been  obtained  by  Chien     and  Feodosyev  -^  . 

In  reference  [6J  Kaplan  and  Flang  present  a  theoretical  and 
experimental  treatment  for  the  clamped  shallow  spherical 
shell  under  uniform  external  pressure.   Their  theoretical 
work  is  based  upon  a  perturbation  solution  of  the  non- 
linear equations,  where  the  non-dimensional  radial  deflec- 
tion at  the  center  of  the  shell  serves  as  the  perturbation 
parameter.   Unfortunately,  solutions  to  the  resulting 
linear  perturbation  equations   could  be  obtained  only  for 
very  shallow  shells.   However,  their  experiments  as  ob- 
served in  reference  [7]   revealed  some  new  and  Interesting 

phenomena  concerning  the  buckling  of  curved  structures. 

[15] 
Archer    ,  using  a  similar  set  of  equations  derived  by 


By  shallow  it  is  implied  that  the  semi-angle  of  opening  -A 
of  the  spherical  segment  is  a  small  angle  j  see  P'ig.  1. 
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Reissner'^"^   ,  extended  and  Improved  the  perturbation 
technique,  and  was  able  to  obtain  solutions  for  a 
greater  range  of  shells  than  Kaplan  and  Fung.   The 
agreement  with  experiments  of  the  buckling  loads  pre- 
dicted from  the  perturbation  solutions  becomes  poorer 
as  the  value  of  the  geometric  shell  parameter  (p  =  -j^) 
increases.   Furthermore,  the  perturbation  solution, 
based  upon  a  one  parameter  expansion,  falls  to  reveal 
the  non-monotone  dependence  of  buckling  loads  on  p  as 
indicated  in  reference  [71.   Simons  ^  proposes  solving 
the  non-linear  equations     by  a  power  series  technique 
similar  to  that  presented  in  this  paper.   His  calcula- 
tions however  cover  only  a  small  range  of  parameters. 

The  physical  problem  that  is  treated  in  this  paper 
is  that  of  a  thin,  shallow  spherical  shell  \inder  the 
action  of  uniform  external  pressure.   This  shell  is 
ass\iraed  fixed  and  clamped  along  the  edge.   Only  rotation- 
ally  symmetric  deformations  will  be  considered.   The 
problem  is  considered  as  one  in  the  non-linear  theory  of 
elasticity  with  small  strains  and  small   but  finite 
deflections.   Non-linear  load  vs  deformation  character- 
istics of  the  shell  is  expected  to  have  the  form  shown 
in  Fig.  2.   Here,   p  represents  the  external  pressure 
and  w  the  radial  deflection  at  some  point,  i.e.  the 
center.   Point  A  represents  the  "upper"  or  "initial" 
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buckling  load.   Point  B  represents  the  "final"  or 
"minimum"  buckling  load. 

The  present  investigation  has  been  successful  in 
obtaining  solutions  including  deflections  and  stresses 
for  a  large  range  of  shell  geometries  for  pressures  from 
zero  to  those  in  the  vicinity  of  point  A  of  Fig.  2. 
Unfortunately  only  meagre  results  have  been  obtained  so 
far  for  the  "post-buckling"  region,  i.e.  the  portion  of 
the  curve  A,  B,  C   in^ig,  2. 

The  teclmique  of  solution  employed  was  that  of 
power  series  expansions  inspired  by  the  success  of 
Friedrichs  and  Stoker    ,  and  Bromberg     in  the  prob- 
lems of  bending  and  buckling  of  circular  plates.   The 
use  of  this  approach  in  the  present  problem,  however, 
was  only  possible  due  to  the  availability  of  a  high- 
speed digital  computer. 
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Section  2.   Derivation  of  Differential  Equations 

The  non-linear  equations  for  the  deformation  of 
spherical  shells  of  s mall  curvature,  based  upon  an 
assumption  of  small  strains  and  linear  stress-strain 

r  o  1 

law,  were  originally  obtained  by  Chien    from  static 
considerations.   These  equations  can  also  be  derived 
from  a  variational  principle.   In  this  case  we  employ 
the  Principle  of  Minimum  Potential  Energy,  and  seek 
Euler's  equations  for  the  functional. 


^1^  TT  =  ffl  Ir'j^i/v  +    jj  pw 


d  S 


where:   TT   is  the  potential  energy;   X'^'^      a-^^Q  the  com- 
ponents of  the  stress  tensor;   e. .   are  the  components 
of  the  strain  tensor;   w   is  the  radial  displacement; 
p  is  the  external  pressure,  positive  inwards.   See 
Fig.  1  for  the  shell  geometry. 

The  functional  (1)  may  be  expressed  in  terms  of 
the  components  of  displacement  through  the  following 
relations. 
Stress  Strain 

,2,  ^11  .  ..  =  _£_ju.  .  ,,.)  ,  ^22  ,  ,..  =  _I^u.  .  veg) 

From  the  assumption  of  rotational  symmetry  all  other  stress 
and  strain  components  vanish, 
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strain  Displacement 

Let   r  be  the  radial  distance  to  any  point  in  the 
shell,  and  R  the  radius  of  the  middle  surface.   The 
exact  non-linear  strain  displacement  relations  for  rota- 
tionally  symmetric  deformations  in  a  spherical  coordinate 

system,  are  then, 

2  2 

u w    1  W5  +  u    ,  ug  -  w 

-   2 

^ -Dv  \       »   _  „   _  u  cot  Q-w  .  1/UcotO-Wx 
(3b)       6^2  =  e^ ^ +  2^ ? ) 

where  subscripts  on  u  and  w  denote  differentiation. 
For  the  case  of  small  finite  deflections,   u  and 


w   (analagous  to  the  Foppl-Karman  non-linear  plate  theory) 

2 

'9 


_     2 
and  shallow  shells  (small  ©,   U5  Is  negligible  compared 


2 
to  Ug  and  wr,   and  u   Is  negligible  compared  to  wg)  , 

equations   (3a, b)   reduce  to 


w-  2 
(i^a)  eg  =  ^(ug  -  w)  +  i(-|~) 


(Ub)  e;==^(^-w) 

(pre 

Equations  {l^a.,h)      are  the  strain  displacement  rela- 
tions that  will  be  employed  in  the  variational  expression. 
It  is  not  clear  to  what  extent  these  equations  are  adeqioate 
for  non-shallow  shells,  without  retaining  additional  terms 
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from  Equations   (3a,b)j  for  example,  see  page  260  of 
reference  3« 

As  is  customary  In  thin  shell  theory  the  following 
assumptions  are  made  concerning  the  distribution  of  the 
displacements  through  the  shell  thickness. 


/  r-S  o    <^   o  o 

(5)  u  =  u-^wg;w  =  w 


where   ^^  =  R-r  and  u   and  w   are  the  displacements 
of  the  middle  surface   ^  =  0. 

Substituting   (2),  (i4.a,b)  and  (5)  into  (1),  and 
performing  the  T,     and  "$  Integrations,  the  potential 
energy  JJ"  may  be  written  as. 


(6) 


rk 


rJ^ 


"  ]  liT^  i  ^"^^' "  1^  "i "  r  "§"§5^  '^' "  \p"^ 


where   e^  and   e-r   are  the  strains  on  the  middle  surface, 
Setting  the  first  variation  of  Equation  (6)  equal  to 
zero,  the  following  equilibrium  equations  result. 


(7a)  o2  =  d_  (5  ^) 

<!>   de    ® 
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^^^e^-i^^^lfs^^-S-e) 


(7b) 

rM-  1  de   de  e  de   ® 


R 


D  1  d^  5  d^  i  d^  g  ^^  _  pQ  ^  0    ^ 


where  o^     and   ov   are  stresses  on  the  middle  surface,  and 


3(l-v^) 


Equations  (7a,b)  are  a  set  of  two  equations  in  the 
three  unknowns,   O^,  a^     and  w.   It  is  also  necessary  to 
satisfy  the  compatibility  equation  for  the  middle  surface, 
This  equation  can  be  obtained  by  eliminating  u  from 
Equations  (l4.a,b)  for   r  =  R.   There  results. 


(8)   i,  eCo-2  -  voS).  o|  .  vo-2  .  I  §  W5  .  I  (^2,"  .  0 


In  addition  to  yielding  the  differential  equations 
(7a,b),  the  variational  problem  leads  to  the  following 
boundary  conditions 


f 


(9)  _A 
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The  first  equation  of  equilibriuin  (7a)  is  identi- 
cally satisfied  if  we  set 

(10a,b)  ^=lhl»   ^|=lh^   • 

Furthermore,  let  us  introduce  the  more  convenient 
variables 


(11)  a  =i^  +  e 

^      '  R  d© 


and 


a  =  Aa      :      Q  =   Ae   :   k  =  5- 

3(l-v^) 

(12)       Y  =  (ic)^/^  4^  Y  ;   P  =  2E(k)l/2(|)2^   . 


=  (k)-l/2_^2  R 


Geometrically,  the  independent  variable   a  defined 
in  Equation  (11)  is  the  angle  that  the  tangent  plane  to 
the  deformed  surface  makes  with  the  horizontal  plane, 
see  Fig.  1. 

The  parameter   p,   hereafter  referred  to  as  the 
"geometrical  parameter",  is  proportional  to  the  ratio 
of  the  shell  rise  and  the  thickness.   The  shell  rise 
is  defined  as  the  height  of  the  center  of  the  shell 
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above  the  horizontal  plane  through  the  edge.   p  Is 
also  roughly  proportional  to  the  ratio  of  membrane  and 
bending  stiffness.   The  relation  between   p  and  Kaplan 
and  Fung's  geometrical  parameter   X   Is, 


(12a)  p  =  (2)"^/^  X^ 


Substituting  (10a, b),  (11)  and  (12)  Into  (7a, b) 
and  (8),  and  Integrating  (7b)  once  (employing  the  fact 

dW  \ 

that  from  rotational  symmetry  -rg  =0  at   ©  =  0;   the 
following  differential  equations  are  obtained. 


(13a, b)  La  =  p(aY  +  P©^)  J   Ly  =  p( e^  -  a^) 


where   L  Is  the  linear  differential  operator. 


(ii^'  ^-^h^h^ 


The  boiondary  conditions  corresponding  to  the  dif- 
ferential equations  (13a,b)"are  as  follows:   At  the 
center,   9=0,   from  rotational  symmetry 

a  =  0   ;    u  =  0 

From  Equations  (2),  (i].a,b)  and  (10a,b)  this  Is  equiv- 
alent to 

(15)  a(0)  =  0   ;    y(0)  =  0 
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At  the  edge,   ©  =  1,   we  assume  that  the  plate  is 
held  fixed  and  clamped.   This  condition  can  also  be 
considered  as  simulating  the  boundary  separating  the 
region  of  dimpling  from  that  of  uniform  contraction 
of  a  complete  sphere.   Therefore,  the  conditions  are 

a  =  l   :   u=0   :   w=0 


This  is  equivalent  to 


(16)  a(l)  =  1   i     Ya  - 


vy|      =  0 
J  9=1 


It  should  be  observed  that  if  conditions  (l5)  and  (16) 
are  imposed  then  the  boundary  conditions  (9),  resulting 
from  the  variational  problem,  are  identically  satisfied. 

Thus,  the  complete  formulation  of  the  problem  con- 
sists of  the  differential  equations  (13a,b)  subject  to 
the  boundary  conditions  (15)  and  (16),   These  equations— 
a  pair  of  second  order,  coupled,  ordinary,  non-linear 
differential  equations  —are  quite  similar  to  the  equa- 
tions governing  the  non-linear  deflections  of  circular 
plates,   *   *    except  for  the  addition  of  extra  terms 
resulting  from  the  initial  curvature.   However  the 
physical  phenomena  that  these  equations  describe  are 
quite  different.   Mathematically,  the  buckling  problem 
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for  the  circular  plate,  or  the  bending  problem  differ 
in  their  parametric  dependence  from  the  spherical 
shell  problem.   The  plate  problems  depend  upon  one 
parameter,  either  a  reduced  edge  load  (the  "thrust 
ratio",  see  Reference  [10]) or  a  reduced  lateral  load, 
see  Reference  [11].   While  on  the  other  hand,  the  equa- 
tions governing  the  shell  problem  depend  upon  two  para- 
meters,  P  and   p.   Prom  this  point  of  view,  it  would 
seem  that  to  obtain  a  more  reasonable  perturbation 
solution,  for  the  shell  problem  a  two  parameter  expan- 
sion should  be  attempted. 
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Section  3«   Technique  of  Solution 

The  technique  of  solution  to  oe  employed  for  the 
boundary  value  problem  posed  by  i)quatlons  (13a ,b),  (l5), 
and  (16)  is  to  expand  the  dependent  variables   a  and 

Y  in  power  series  in  6.   This  approach  has  been  used 

l9] 
in  non-linear  circular  plate  problems  by  Way   , 

Priedrichs  and  Stoker     and  Bromberg 

Prom  (13a,b)  and  the  boundary  conditions  at   9=0 

the  power  series  for   a  and   y  roust  contain  odd  powers 

only   so   that 


00  o       T  CO  o       , 

(17)  a  =  i:V         '    r  =  rY  e^""^ 

n=l  n=l 


By  substitution  of  these  power  series  into  (13a,b)  the 

following  recurrence  relations  for  the  coefficients  a 

n 

and      Y       are   obtained, 
n 


a.^  =   p/8(a^Yi  +  P)    ;      Y2  =   P/8(l   -  a^)    ; 
(18)        /       for     n  >  2 

^n  =  i|n(S-l)^5^i^n-i    '      ^n  =   "   i^(g-l)^|l^i^n-l 

Again,  by  insertion  of  ;;he  power  series  into  Equation 
(16)  the  boxmdary  conditions  at   9  =  1   can  be  written  as 
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00  00 


3"  a  =  1   ;  ri  (2n-l-v)Y„  =  0 
n=l  n-1 

Since  it  is  possible,  through  the  recurrence  relations 
(18),  to  express  all   a^  and   y^   i^  terms  of  a^  and 
Y,   the  above  conditions  constitute  a  pair  of  simul- 
taneous equations  in  a^   and   y^^   i.e. 

OD  CO 

(19)  F{a     r   )    zL!v  1  =  0  J  G(a^,Yi)  =  UlSn-Cl+v)  ]  y^=  < 
n=l  n=l 

The  problem  thus  reduces,  for  given   P  and   p,   to  solv- 
ing the  two  simultaneous  equations  (19)  for  the  unknowns 
a,   and   y-i  »   Once  the  root  or  roots  of  these  equations 

are  obtained  all  the   a   and   y   and  hence  all  the 

n        n 

stresses  and  displacements  are  then  determined;  (see 
Appendix  A) . 
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Section  k»      Computational  Technique 

The  computational  problem  is,  for  a  given  set  of 
parameter  values   (P,p)   to  determine  the  roots   ^*^i»Y]^^ 
satisfying  (19).   For  computing  purposes,  it  is  necessary 
to  replace  the  exact  equations  (19)  by  the  approximate 
conditions 

N  N 

(19')  |F(a^,Yi)l^lZ  V  ^I<1  '  |G(a^,Yi)|*IZ  t2n-(l+v)]Y„|< 

where  'V      is  a  small  niimber  and   N   is  an  integer  deter- 
mined on  the  basis  of  a  convergence  criterion  for  the 
series.   In  the  majority  of  the  computations  '^      was 
chosen  as  1  x  lO"^.   N   is  never  greater  than  99. 

For  a  fixed   P  and   p  the  procedure  starts  from 
a  reasonable  guess  for   (a^,Y^),   say   ^'^t*Y-i)  •   Then 
from  (18)   a   and  y       corresponding  to   (a^,Y-)  ) 
are  computed  for  n  <  N,  As  a  convergence  criterion 
for  the  series  (17),  it  is  prescribed  that  the  sum  of 
the  absolute  magnitudes  of  the  last  three  computed  a 

and   Y   be  less  than  a  certain  small  number,   e.   That 
'n  ' 

is,  the  computation  of  a   and  y„   in  (16)  ceases  for 
n  =  N   such  that 


(20)  I^mI  +  I^n.iI  +  Ivz'^^  '  1%I"*'1yn-i1''I%-2I<^ 
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In  the  majority  of  the  computations   e  was  set  equal  to 
1  X  10"^. 

When  the   a   and   y   satisfy  Equation  {20),    they 
are  inserted  into  the  expressions  for   F(a^,Y-,)   and 
G{a   ,y).      If  F  and   G   satisfy  Equation  (19 ') ,  the 
stress  and  displacement  distribution  throughout  the  shell 
are  computed  from  the  relations  listed  in  Appendix  A, 
In  general,  however,   P  and   G  do  not  satisfy  Equation 
(19' )»  and  it  is  necessary  to  obtain  an  improved  approx- 
imation to   (a^,Yi).   This  is  done  systematically  by 
means  of  the  Newton-Raphson  iteration  process.   The 
partial  derivatives  of  F  and  G  in  the  Newton- 
Raphson  process  are  approximated  by  forward  difference 
quotients.   This  iteration  process  is  repeated  iintil 
Equation  (19')  is  satisfied.   The  displacements  and 
stresses  are  then  computed  from  the  equations  of  Appen- 
dix A.   The  number  of  iterations  required  for  convergence 
varies  with  the   (P,p)   selection  and  the  Initial  guess 
(cL.,Y-|)  .   For  most  of  the  cases  encountered,  two  to 
six  Newton-Raphson  iterations  were  necessary.   The  num- 
ber of  coefficients,   N,   that  were  required  for  con- 
vergence varied  with  the   (P,p)   selection;  Table  I 
gives  some  representative  values.   The  computation  of 
coefficients,  the  iteration  process,  and  the  calcula- 
tion of  stresses  and  displacements  are  all  performed 
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on  the  Univac  at  New  York  University. 

The  method  for  calculating  with  fixed   P  and   p 
has  thus  been  described.   The  computations  proceed,  for 
fixed   p,   by  obtaining  a  sequence  of  roots   ^^i*^^) 
for  a  corresponding  sequence  of  increasing  values  of 
P.   A  first  guess  for   (a^,Y^)   is  determined  by  selec- 
ting a  low  value  for   P,   usually   P  =  .1.   For  such 
low  values  of   P  the  deformed  state,  and  hence  the 
first   (a-,Y-,)   guess,  should  be  nearly  spherical.   The 
Newton-Raphson  process  then  rapidly  yields  the  root 
(a  ,Yn )   to  the  required  accuracy.   The  initial  guess 
of  the  roots  for  successivly  higher  values  of   P   is 
obtained  by  extrapolation  of  the  preceding  roots  in 
the  sequence.   This  is  done  with  the  aid  of  a  plot  of 
a-   vs   Y-]   as  shown  in  Fig.  3«   This  procedure  is 
terminated  for  either  of  two  reasons:   the  Newton- 
Raphson  iteration  process  fails  to  converge  to  a  root 
in  a  large  number  of  iterations  (usually  10);  or  the 
coefficients  fail  to  converge  according  to  Equation 
(20).   As  an  example  of  the  first  reason  for  termi- 
nation, with   p  =  32   and   P  =  1.25   the  procedure  con- 
verged in  3  iterations.   However,  with  P  =  1.3   there 
was  no  indication  of  convergence  in  10  Newton-Raphson 
1 


The  authors  are  indebted  to  Miss  M.  Weissner  and  Messrs. 
A.  Paster  and  S.  B.  Weinberg  for  their  aid  in  coding  and 
programming  this  problem. 
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iterations.   An  Illustration  of  the  second  reason  for 
termination  is  for   p  =  26.1   and   P  =  1.55.   Here  con- 
vergence occured  in  2  iterations  and  1^.0   terras  of  the 
power  series.   But,  for   P  =  1.6,  99   terms  were  computed 
with  Equation  (20)  not  being  satisfied.   In  fact,  the 
terms  were  steadily  growing  in  absolute  magnitude. 

The  criterion  for  the  initial  buckling  that  is 
used  in  this  investigation  is  based  on  this  terminating 

value  of   P,   Direct  application  of  the  "classical" 

dP 
criterion  (the  load  at  which  "^r:  =  0)  ^or   buckling 

seems  Impracticable  with  the  present  technique  of 
solution.   Instead,  the  niamerlcal  equivalent  of  the 
buckling  phenomena  is  used  as  the  criterion  for  de- 
termining the  buckling  load,  i.e.  the  failure  of  the 
computing  process  to  produce  a  neighboring  solution 
corresponding  to  a  small  incrameht  in  the  load.   If 
P,   represents  the  last  value  of   P  for  which  con- 
vergence occurs,  and   Pp  represents  the  next  value  of 
P  attempted,  then  the  initial  buckling  load,   P^_» 

is  assumed  to  lie  in  the  Interval   P,  <  P   <  P^, 

1  —  cr  —  2 
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Section  5»   Discussion  of  Results 

Before  proceeding  to  a  discussion  of  the  results  of 
our  calculations,  a  brief  review  of  the  experimental  re- 
sults of  Kaplan  and  Pung  and  the  interpretation  given 
in  reference  [7]  will  be  presented. 

Kaplan  and  Fung  constructed  a  series  of  shallow  spher- 
ical domes  of  varying  radii  and  thickness.   These  domes  were 
placed  in  a  clamping  jig  and  loaded  under  uniform  external 
pressure.   In  addition  to  m.easuring  the  initial  and  final 
buckling  pressures,  they  measured  the  vertical  deflections 
at  various  points  for  several  different  values  of  the 
pressure . 

For  low  values  of  the  geometric  parameter   p   (or 
equlvalently  Kaplan  and  Fung's  X     as  related  by  Equation 
12a)  the  experiments  showed  that  the  form  of  the  deflec- 
tions are  such  that  the  maximum  deflection  is  at  the 
center  of  the  shell,  and  decreases  steadily  toward  the 
edge.   This  will  be  called  mode  of  deformation   I,   As 
p  increases,  it  was  observed  that  the  peak  deflection 
at  the  center  gradually  flattened  out,  until  at   p 
equal  approximately  21,  the  maximum  deflection  occurs 
at  some  point  between  the  center  and  edge;  so  that 
there  is  also  a  local  minimum  point  at  the  center. 
This  will  be  called  mode  of  deformation  II,   With  a 
further  increase  in   p  the  maximum  point  moved 
toward  the  edge  until  approximately   p  =  SS*      Then  the 
shell  deformed  in  a  third  mode.   This  mode  has  two 
local  maximum  points,  one  at  the  center  and  the  other 
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between  the  center  and  the  edge,  see  Pig.  l6  of  Refer- 
ence 6. 

When  the  experimental  buckling  loads  are  plotted 

[7] 

against   p  rather  than  X,   it  is  found    that  there 
is  an  Interesting  peaking  behavior  in  these  loads  as  a 
function  of   p,   see  P±g,   I4..   In  Fig.  l^   the  crosses  and 
circles  indicate  the  experimental  values  of  Kaplan  and 
Fung,  while  the  solid  line  indicates  the  apparent  trend 
of  the  experimental  data.   It  is  significant  to  note 
that  these  peaks  in  the  buckling  loads  occur  at  those 
values  of   p  at  which  a  change  of  mode  of  deformation 
had  been  observed. 

In  proceeding  with  our  numerical  solutions  on  the 
Univac  we  have  in  a  real  sense  conducted  "experiments" 
on  shallow  spherical  shells.   For  each  given  value  of 
p   (i.e.  a  particular  shell  configuration)  deflection 
and  stress  "measurements"  are  made  for  increasing  val- 
ues of  load   P.   As  in  a  physical  experiment,  the  ini- 
tial buckling  load  is  considered  to  be  reached,  when 
for  a  small  increment  in   P  no  neighboring  solution 
to  the  previously  attained  solution  can  be  found,  i.e. 
the  iteration  process  suddenly  fails  to  converge. 

For  low  values  of   p,   p  <  23,  the  numerical  re- 
sults confirm  Kaplan  and  Fung's  experiments.   In 
Figures  5a,b,c,d  are  plotted  some  values  of  the  deflec- 
tions and  stresses  as  a  function  of  9.   The  plots  are 
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for  increasing  load   P  at  a  constant  value  of   p.   In 
this  case   p  =  11.312  was  selected  as  representative 
of  the  mode  of  deformation  of  type,  or  region  I.   In 
Fig.  5a  the  radial  deflection  w  can  be  seen  to  peak 
at  the  center  and  fall  off  towards  the  edge.   The  dotted 
lines  indicating  the  results  from  the  perturbation 
solution  of  Kaplan  and  Fung*-   .   It  can  also  be  noted 
from  Figures  5b,c,d  that  the  stresses  exhibit  a  similar 
behavior. 

For  increasing  values  of   p,   the  graphs  of  these 
functions  tend  to  flatten  out  in  the  central  region  of 
the  shell   (8  =  0).   This  process  continues  xintil 
p  =  23.5,   see  Fig.  6  where  only  the  results  for  the 
radial  deflection  are  presented.   For  this  value  of   p 
and  for  low  values  of  load,  say   P  <  .Ij.,   the  shell 
deforms  in  mode  I  with  maximiim  deflection  at  the  center. 
However,  for  higher  load  the  mode  of  deformation  grad- 
ually changes  into  type  II,  with  the  msiximxim  deflection 
located  somewhere  between  the  center  and  the  edge  (in 
this  case   9Jy.3).   It  is  to  be  noted  from  Fig.  i|  that 
the  shell   p  =  23.5  falls  in  the  peaking  region  of  the 
buckling  loads.   Thus,  the  mode  of  deformation  for  con- 
stant  p  may  change  with  increasing  load.   This  clearly 
Illustrates  the  two  parameter  dependence  of  the  charac- 
ter of  the  solutions. 
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In  Figiires  7a,b,c,d  are  shown  plots  of  deflections 
and  stresses  for  a  typical  shell    p  =  32,   of  mode  II 
of  deformation  for  increasing  values  of  load.   Notice, 
that  the  longitudinal  membrane  stress  Fig.  7b  also  has 
a  maximum  located  between  the  center  and  edge  (compare 
to  Pig.  5b).   The  behavior  of  the  bending  stresses  for 
mode  II  (see  Figures  7c, d),  shows  a  marked  difference 
from  that  of  mode  I  (see  Figures  5c,d).   For  mode  I, 
the  bending  stresses  fall  off  in  a  smooth  manner  from 
the  center  to  the  edge,   while  for  mode  II  there  is  a 
local  minimum  near   0  =  ,6  and  then  the  stresses  change 
rapidly  towards  the  edge.   This  is  reminiscent  of  the 
boundary  layer  phenomenon  foxind  by  Friedrichs  and  Stoker 
for  the  buckling  of  circular  plates. 

A  further  increase  in   p  results  in  a  flattening 
of  the  deflection  and  stress  distributions  in  the  centrail 
portion  of  the  shell,  and  a  movement  of  the  maximum 
point  towards  the  edge.   Mode  II  persists  until   p^52 
when  a  transition  case  similar  to   p  =  23.5   is  detected. 
However  this  time  the  transition  is  from  mode  II  to  mode 
III.   In  the  radial  deflection  plot,  the  mode  III  form 
is  barely  visible.   However  it  can  clearly  be  observed 
in  the  bending  stress  graphs.   It  is  not  until  large   p, 
say   p  =  176,   that  there  is  a  sharp  appearance  of  mode 
III  in  the  radial  deflection  graph,  (see  Fig.  8).  Again 
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from  Fig.  li   it  can  be  observed  that  the  transition  of 
modes  occurs  at  the  peaking  of  the  buckling  loads. 

In  Figures  9a,b,c,d  the  deflections  and  stresses 
are  again  plotted  against   6,   but  this  time  for  con- 
stant  P  and  increasing   p.   In  Fig.  9a  the  motion  of 
the  maximum  deflection  point  towards  the  edge  with  in- 
creasing  p  is  clearly  indicated.   This  is  again 
reminiscent  of  boundary  layer  phenomena.   In  the  plots 
for  the  longitudinal  membrane  stress  (Fig.  9b)  and  in 
particular  the  bending  stresses  (Fig.  9c,d)the  change 
in  mode  shape  is  clearly  illustrated. 

The  load  vs  deflection  at  the  center  character- 
istics for  various  values  of  p  are  plotted  in  Fig.  10 
and  compared  with  the  theoretical  and  experimental 
results  obtained  by  Kaplan  and  Fung.   For  low  values 
of   p  the  agreement  with  the  experiments  is  only  fair, 
while  for  increasing   p  the  agreement  between  the 
physical  gmd  numerical  experiments  becomes  uniformly 
better.   This  is  probably  caused  by  the  fact  that  for 
shallow  shells  (small   p)   the  effect  of  initial  imper- 
fections are  more  pronounced,  and  that  the  shallower 
the  shell  the  more  difficult  it  is  to  properly  measure 
the  radius. 

Initial  Buckling  Loads:   The  results  of  the  numer- 
ical coiiiputation  of   P    as  well  as  the  experimental 
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results  of  Kaplan  and  Fiing  are  presented  in  Fig.  l\.» 

The  intervals,   PtPo,   in  which  the   P    lie  are 
'12'  cr 

indicated  by  short  vertical  lines.   It  can  be  observed 

that  the  computed   P    follow  the  same  general  trend 

as  the  experimental  data,  indicating  a  peaking  behavior 

as  a  function  of   p.   The  agreement  between  experiment 

and  theory  is  least  satisfactory  in  the  region  of 

peaking  [18  <  p  <  25,  50  <  p  <  6o].  In  this  region  a 

small  error  in  the  experimental  determination  of   p 

may  lead  to  a  large  error  in   P    for  a  particular 

cr 

shell.   It  can  also  be  observed  from  Fig.  I|.  that  in 
some  cases  the  numerical  results  predict  greater   P 
than  the  experiments,  while  in  other  cases  it  is  just 
the  opposite.   There  is  no  apparent  trend  In  this 
respect. 
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Section  6.   Post  Buckling  Behavior 

The  post-buckling  region  Is  defined  as  the  region 
A,B,C,  of  Fig.  2,  i.e.  the  behavior  of  the  shell  after 
P  =  P   .   The  primary  purpose  of  investigating  this 
region  is  to  determine  the  lowest  value  of   P,(P  =  P^) 
for  which  there  exist  non-ixnique  solutions  to  the 
boundary  value  problem  described  by  Equations  (13a,b^ 
(l5)  and  (16).   The  procedure  previously  outlined  does 
not  enable  us  to  continuously  obtain  solutions  beyond 
point  A.   However,  starting  from  a  solution  in  the 
neighborhood  of  point   C  we  may  apply  the  previous 
procedure  with  a  decreasing  sequence  of   P  values  to 
obtain  the  point  B. 

To  obtain  a  solution  in  the  neighborhood  of  point 
C  we  select  a  sufficiently  large  value  of   P   (say 
P  =  2.2).   For  a  "flat"  shell  the  load  vs  deformation 
characteristics  must  be  monotone;  I.e.  there  Is  no 
buckling.  .Therefore  for  a  sufficiently  small  value 
of  p  it  is  possible  to  obtain  a  solution,  using  the 
previous  procedure,  for  this  value  of  P.   Then  start- 
ing with  the  solution  for  this  value  of   p  we  obtain 
a  sequence  of  solutions  for  an  increasing  sequence  of 
p  values  and  fixed  P  (see  Fig,  11,  p  =  .8,  i^-.8,  S.ij., 
7.2i|.).   This  procedure  is  analogous  to  the  original 
procedure  where  the  roles  of   P  and   p  are  interchanged. 
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We  note  how  for  increasing   p  values  the  character- 
istics depart  more  and  more  from  the  initial  linear 
relation  and  approach  the  buckling  type  curve. 

A  solution  is  thus  obtained  for   P  =  2.2  and 
p  =  11.312.   This  load  is  higher  than  that  previously- 
used  and  corresponds  to  much  larger  deflections. 
Keeping   p  fixed  a  decreasing  sequence  of   P  values 
is  selected  in  an  attempt  to  locate  point   B.   Un- 
fortunately, convergence  of  the  computations  becomes 
increasingly  poor  as   P  decreases.   Solutions  are 
obtained  only  for  the  range   1.8  <  P  <  2»2,      see  Fig. 
(11).   The  initial  portion  of  the  curve  for   p  =  11.312 
is  also  shown  in  Fig,  11.   Fig.  12  represents  the 
radial  deflection  over  the  shell  for  one  value  of   P 
in  the  post-buckling  region. 
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Section  7»   Conclusions 

Numerical  solutions  of  the  non-linear  equations 
governing  the  deflections  of  clamped  shallow,  spheri- 
cal shells  under  a  uniforru  external  pressure  have 
been  obtained  by  power  series  technique.   The  compu- 
tations were  conducted  on  the  Univac  at  N.  Y.  U, 
These  solutions  yield  stress  and  deflection  distribu- 
tions over  the  shell  as  well  as  estimates  of  the 
Initial  buckling  loads  for  a  wide  range  of  values  of 
the  shell  parameters.   The  two  parameter  dependence 
of  the  solutions  and  the  various  modes  of  deformation 
have  been  clearly  shown.   Only  a  limited  number  of 
solutions  for  the  post-buckling  region  have  been  ob- 
tained.  It  is  anticipated  that  some  Improvements  in 
the  numerical  method  may  yield  additional  results  in 
this  region. 

The  authors  wish  to  express  their  appreciation 
to  Profesgor  J.  J.  Stoker  for  suggesting  this 
investigation. 
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APPENDIX  A 
Once  the  roots   ^^i*T-,  )   are  obtained  for  a  given 
set  of  parameter  values   (P,p),   the  stresses  and  de- 
flections may  be  computed  from  the  following  relations: 
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